#---- Replication script 2/4----

#Turnbull-Dugarte
#"Liberal in the sheets [...]"

###This script produces all analysis from Study II

remove(list=ls())

#install.packages("pacman")
library(pacman)

p_load(here, tidyverse, dplyr, haven, cregg, modelsummary, rcartcolor, jtools, broom, sdamr,
       distributional, list, ggpubr, sruvey,  ggdist, patchwork, estimatr, margins, grf, starbility,
       miceadds, modelr)


col5 <- c("navy","navy")
col4 <- c("purple", "purple4",
          "darkorange", "darkorange4")
col2 <- c("#FC4C02", "purple")
col2b <- c("purple","#FC4C02")

###UK CLEAN UP##

tinder_conjoint_UK <- readRDS("data/study2_UK.rds")
tinder_conjoint_UK <-tinder_conjoint_UK%>% 
  mutate(BI= as.factor(BISEXUAL),
         TRANS= as.factor(TRANS),
         CONGRUENCE= as.factor(CONGRUENCE),
         FLAG= as.factor(NATIONALITY),
         DISTANCE= as.factor(DISTANCE),
         RELATIONSHIP= as.factor(RELATIONSHIP),
         PROS= as.factor(PROS),
         resp_gender=as.factor(gender),
         PARTY=as.factor(PARTY),
         experimental_group=as.factor(experimental_group),
         EDUCATION_SLIM= as.factor(EDUCATION_SLIM),
         JOB_SLIM= as.factor(JOB_SLIM),
         ETHNICITY= as.factor(ETHNICITY),
         EMOTION= as.factor(EMOTION))

tinder_conjoint_UK$BI<- factor(tinder_conjoint_UK$BI,
                               levels = c("FALSE","TRUE"),
                               labels = c("Heterosexual\nprofile", "Bisexual\nprofile"))

tinder_conjoint_UK$TRANS<- factor(tinder_conjoint_UK$TRANS,
                                  levels = c("FALSE","TRUE"),
                                  labels = c("Cis.", "Trans."))

tinder_conjoint_UK$FLAG<- factor(tinder_conjoint_UK$FLAG,
                                 levels = c("No flag", "UK flag", "UK-EU flag", "English flag", "English-EU flag", "EU flag"))


tinder_conjoint_UK$derecha<- factor(tinder_conjoint_UK$derecha_num,
                                    levels = c(0, 1),
                                    labels = c("Left-wing", "Right-wing"))

tinder_conjoint_UK <- tinder_conjoint_UK |> 
  mutate(derecha_mujer= factor(case_when(derecha== "Right-wing" & gender=="a woman" ~ "Right-wing woman",
                                         derecha== "Right-wing" & gender=="a man" ~ "Right-wing man",
                                         derecha== "Left-wing" & gender=="a woman" ~ "Left-wing woman",
                                         derecha== "Left-wing" & gender=="a man" ~ "Left-wing man",
  )))

tinder_conjoint_UK <- tinder_conjoint_UK |> 
  mutate(queer= factor(case_when(sex_id== "0" ~ "Heterosexual",
                                 sex_id== "1" | sex_id== "2" |sex_id== "3" |sex_id== "4" ~ "LGBTQ+"
  )))

tinder_conjoint_UK <- tinder_conjoint_UK |> 
  mutate(FARRIGHT= factor(case_when(PARTY== "#VoteReformUK" ~ "Far-right supporter",
                                    PARTY== "#VoteLibDem" | PARTY== "#VoteLabour" |PARTY== "#VoteConservative" |PARTY== "#VoteGreens" ~ "Not far-right supporter"
  )))


tinder_conjoint_UK$PARTY<- factor(tinder_conjoint_UK$PARTY,
                                  levels = c("#VoteReformUK","#VoteLibDem", "#VoteLabour", "#VoteConservative", "#VoteGreens"),
                                  labels = c("Reform UK", "Lib Dem", "Labour", "Conservatives", "Greens"))
tinder_conjoint_UK$PARTY <- relevel(tinder_conjoint_UK$PARTY, ref = "Conservatives")


table(tinder_conjoint_UK$PARTY, tinder_conjoint_UK$FARRIGHT)
table(tinder_conjoint_UK$selected)

## SPAIN CLEAN UP###
tinder_conjoint_ES <- readRDS("data/study2_ES.rds")
tinder_conjoint_ES$degree <- ifelse(tinder_conjoint_ES$educ >= 4, "Yes", "No")


tinder_conjoint_ES <-tinder_conjoint_ES%>% 
  mutate(BI= as.factor(BISEXUAL),
         TRANS= as.factor(TRANS),
         CONGRUENCE= as.factor(CONGRUENCE),
         FLAG= as.factor(NATIONALITY),
         DISTANCE= as.factor(DISTANCE),
         RELATIONSHIP= as.factor(RELATIONSHIP),
         VEG= as.factor(VEG),
         resp_gender=as.factor(gender),
         PARTY=as.factor(PARTY),
         experimental_group=as.factor(experimental_group),
         EDUCATION_SLIM= as.factor(EDUCATION_SLIM),
         JOB_SLIM= as.factor(JOB_SLIM),
         ETHNICITY= as.factor(ETHNICITY),
         EMOTION= as.factor(EMOTION))


tinder_conjoint_ES
tinder_conjoint_ES$BI<- factor(tinder_conjoint_ES$BI,
                               levels = c("FALSE","TRUE"),
                               labels = c("Heterosexual\nprofile", "Bisexual\nprofile"))

tinder_conjoint_ES$TRANS<- factor(tinder_conjoint_ES$TRANS,
                                  levels = c("FALSE","TRUE"),
                                  labels = c("Cis.", "Trans."))



tinder_conjoint_ES$derecha<- factor(tinder_conjoint_ES$derecha_num,
                                    levels = c(0, 1),
                                    labels = c("Left-wing", "Right-wing"))

table(tinder_conjoint_ES$derecha, tinder_conjoint_ES$party_id)
tinder_conjoint_ES <- tinder_conjoint_ES |> 
  mutate(derecha_mujer= factor(case_when(derecha== "Right-wing" & gender=="a woman" ~ "Right-wing woman",
                                         derecha== "Right-wing" & gender=="a man" ~ "Right-wing man",
                                         derecha== "Left-wing" & gender=="a woman" ~ "Left-wing woman",
                                         derecha== "Left-wing" & gender=="a man" ~ "Left-wing man",
  )))

tinder_conjoint_ES <- tinder_conjoint_ES |> 
  mutate(queer= factor(case_when(sex_id== "0" ~ "Heterosexual",
                                 sex_id== "1" | sex_id== "2" |sex_id== "3" |sex_id== "4" ~ "LGBTQ+"
  )))

tinder_conjoint_ES <- tinder_conjoint_ES |> 
  mutate(FARRIGHT= factor(case_when(PARTY== "#VotaVox" ~ "Far-right supporter",
                                    PARTY== "#VotaPSOE" | PARTY== "#VotaSumar" |PARTY== "#VotaPP" ~ "Not far-right supporter"
  )))

table(tinder_conjoint_ES$PARTY)
tinder_conjoint_ES$PARTY<- factor(tinder_conjoint_ES$PARTY,
                                  levels = c("#VotaPP","#VotaPSOE", "#VotaSumar", "#VotaVox"),
                                  labels = c("PP", "PSOE", "Sumar", "VOX"))
tinder_conjoint_ES$PARTY <- relevel(tinder_conjoint_ES$PARTY, ref = "PP")



tinder_conjoint_ES <- select(tinder_conjoint_ES, -relationship_type)
tinder_conjoint_UK <- select(tinder_conjoint_UK, -relationship_type)
tinder_conjoint_ES <- mutate(tinder_conjoint_ES, country = "ES")
tinder_conjoint_UK <- mutate(tinder_conjoint_UK, country = "UK")
                                                          
df_full <- bind_rows(tinder_conjoint_ES, tinder_conjoint_UK)
df_hetero <- subset(df_full, queer== "Heterosexual")
table(df_hetero$selected)
table(df_hetero$BI)
df_hetwomen <- subset(df_hetero, resp_gender== "a woman")
table(df_hetwomen$selected)
df_hetmen <- subset(df_hetero, resp_gender== "a man")
table(df_hetmen$selected)

###OVERALL PLOTS###

mm.1 <- mm(df_hetero, selected ~ FLAG + BI + TRANS + CONGRUENCE + ETHNICITY + EMOTION + RELATIONSHIP + DISTANCE + EDUCATION_SLIM + JOB_SLIM  + experimental_group,
             id = ~ ResponseId, h0 = 0.5)

mm.1$lower90 <- mm.1$estimate - (1.645 * mm.1$std.error)
mm.1$higher90 <- mm.1$estimate + (1.645 * mm.1$std.error)
mm.1$lower95 <- mm.1$estimate - (1.96 * mm.1$std.error)
mm.1$higher95 <- mm.1$estimate + (1.96 * mm.1$std.error)
mm.1$lower99 <- mm.1$estimate - (2.58 * mm.1$std.error)
mm.1$higher99 <- mm.1$estimate + (2.58 * mm.1$std.error)

mm.1<-subset(mm.1, feature %in% c("BI")) 

mmplot1<- mm.1 %>% 
  ggplot(aes(level, color=level)) +
  scale_fill_manual(values = col5) +
  scale_color_manual(values = col5)+
  geom_segment(aes(x=lower90, xend=higher90, y=level, yend=level), size=4, alpha=1)+ 
  geom_segment(aes(x=lower95, xend=higher95, y=level, yend=level), size=4, alpha=.6)+ 
  geom_segment(aes(x=lower99, xend=higher99, y=level, yend=level), size=4, alpha=.3)+
  facet_grid(feature~., scales='free', space='free') +  
  geom_point(aes(x=estimate, y=level), fill="white", color="grey67", size=2.5,  pch=21)+
  labs(y=NULL, subtitle='Marginal means (MM)', x="")+
  theme_ggdist()+
  theme(legend.position = "none") +
  xlim(.35,.5)+
  geom_label(
    aes(label = sprintf("%0.2f", round(estimate, digits = 2)), y = level, x = estimate),
    nudge_y = .16, size=2.5, fontface = "bold")+
  theme(strip.background = element_blank(),
        strip.text = element_blank(),
        plot.subtitle = element_text(face="bold", hjust=.5)) + 
  geom_vline(xintercept = 0.44, colour="gray", linetype="dashed") ##.46 median acceptance


amce.1 <- cj(df_hetero, selected ~ FLAG + BI + TRANS + CONGRUENCE + ETHNICITY + EMOTION + RELATIONSHIP + DISTANCE + EDUCATION_SLIM + JOB_SLIM  + experimental_group,
               estimate="amce", id = ~ ResponseId)
amce.1$lower90 <- amce.1$estimate - (1.645 * amce.1$std.error)
amce.1$higher90 <- amce.1$estimate + (1.645 * amce.1$std.error)
amce.1$lower95 <- amce.1$estimate - (1.96 * amce.1$std.error)
amce.1$higher95 <- amce.1$estimate + (1.96 * amce.1$std.error)
amce.1$lower99 <- amce.1$estimate - (2.58 * amce.1$std.error)
amce.1$higher99 <- amce.1$estimate + (2.58 * amce.1$std.error)

amce.1<-subset(amce.1, feature %in% c("BI")) 
amce.1<-subset(amce.1, level %in% c("Bisexual\nprofile")) 

amceplot1<- amce.1 %>% 
  ggplot(aes(level, color=level)) +
  scale_fill_manual(values = col5) +
  scale_color_manual(values = col5)+
  geom_segment(aes(x=lower90, xend=higher90, y=level, yend=level), size=4, alpha=1)+ 
  geom_segment(aes(x=lower95, xend=higher95, y=level, yend=level), size=4, alpha=.6)+ 
  geom_segment(aes(x=lower99, xend=higher99, y=level, yend=level), size=4, alpha=.3)+
  facet_grid(feature~., scales='free', space='free') +  
  geom_point(aes(x=estimate, y=level), fill="white", color="grey67", size=2.5,  pch=21)+
  labs(y=NULL, subtitle='Average marginal component effect (AMCE)', x="")+
  theme_ggdist()+
  theme(legend.position = "none") +
  xlim(-.12,.02)+
  geom_label(
    aes(label = sprintf("%0.2f", round(estimate, digits = 2)), y = level, x = estimate),
    nudge_y = .09, size=2.5, fontface = "bold")+
  theme(strip.background = element_blank(),
        strip.text = element_blank(),
        plot.subtitle = element_text(face="bold", hjust=.5),
        axis.text.y = element_blank(),
        axis.ticks.y = element_blank())+ 
  geom_vline(xintercept = 0, colour="gray", linetype="dashed") ##.46 median acceptance

mmplot1+amceplot1+
plot_annotation(caption = 'Confidence intervals at 99%, 95% & 90% (N:31,240)',
                title = 'Study 2: Results from visual conjoint experiment') & 
  theme(plot.title = element_text(face="bold"))
ggsave(filename = here("figures/Figure7.png"), height=3.5, width=7.5, dpi = 500)
ggsave(filename = here("figures/tifs/Figure7.tif"), height=3.5, width=7.5, dpi = 500)


byGEN.mm<- cj(df_hetero, selected ~ FLAG + BI + TRANS + CONGRUENCE + ETHNICITY + EMOTION +RELATIONSHIP + DISTANCE + EDUCATION_SLIM + JOB_SLIM  + experimental_group,
                id = ~ ResponseId, estimate="mm",by = ~ resp_gender)

byGEN.mm$lower90 <- byGEN.mm$estimate - (1.645 * byGEN.mm$std.error)
byGEN.mm$higher90 <- byGEN.mm$estimate + (1.645 * byGEN.mm$std.error)
byGEN.mm$lower95 <- byGEN.mm$estimate - (1.96 * byGEN.mm$std.error)
byGEN.mm$higher95 <- byGEN.mm$estimate + (1.96 * byGEN.mm$std.error)
byGEN.mm$lower99 <- byGEN.mm$estimate - (2.58 * byGEN.mm$std.error)
byGEN.mm$higher99 <- byGEN.mm$estimate + (2.58 * byGEN.mm$std.error)
byGEN.mm$lower999 <- byGEN.mm$estimate - (3.29 * byGEN.mm$std.error)
byGEN.mm$higher999 <- byGEN.mm$estimate + (3.29 * byGEN.mm$std.error)

byGEN.mm<-subset(byGEN.mm, feature %in% c("BI")) 

mmplot2<- ggplot(byGEN.mm, aes(y = level, x = estimate, color = BY)) +
  scale_shape_manual(values = c(10, 18)) +
  scale_fill_manual(values = col2b) +
  scale_color_manual(values = col2b) +
  geom_vline(xintercept = 0.35, colour = "purple", linetype = "dashed")+
  geom_vline(xintercept = 0.56, colour = "#FC4C02", linetype = "dashed")+
  geom_segment(data = subset(byGEN.mm, BY == "a woman"),
               aes(x = lower90, xend = higher90, y = level, yend = level),
               size = 4, alpha = 1,  position = position_nudge(y = .2)) + 
  geom_segment(data = subset(byGEN.mm, BY == "a woman"),
               aes(x = lower95, xend = higher95, y = level, yend = level),
               size = 4, alpha = .6,  position = position_nudge(y = .2)) + 
  geom_segment(data = subset(byGEN.mm, BY == "a woman"),
               aes(x = lower99, xend = higher99, y = level, yend = level),
               size = 4, alpha = .3,  position = position_nudge(y = .2)) +
  geom_point(data = subset(byGEN.mm, BY == "a woman"),
             aes(x = estimate, y = level), fill = "white", color = "grey67", size = 2.5, pch = 21,  
             position = position_nudge(y = .2))+
  geom_segment(data = subset(byGEN.mm, BY == "a man"),
               aes(x = lower90, xend = higher90, y = level, yend = level),
               size = 4, alpha = 1) + 
  geom_segment(data = subset(byGEN.mm, BY == "a man"),
               aes(x = lower95, xend = higher95, y = level, yend = level),
               size = 4, alpha = .6) + 
  geom_segment(data = subset(byGEN.mm, BY == "a man"),
               aes(x = lower99, xend = higher99, y = level, yend = level),
               size = 4, alpha = .3) +
  geom_point(data = subset(byGEN.mm, BY == "a man"),
             aes(x = estimate, y = level), fill = "white", color = "grey67", size = 2.5, pch = 21) +
  geom_label(data = subset(byGEN.mm, BY == "a woman"),
             aes(label = sprintf("%0.2f", round(estimate, digits = 2)), y = level, x = estimate),
             position = position_nudge(y = .05), size = 2, fontface = "bold"
  ) +
  geom_label(data = subset(byGEN.mm, BY == "a man"),
             aes(label = sprintf("%0.2f", round(estimate, digits = 2)), y = level, x = estimate),
             position = position_nudge(y = -.15), size = 2, fontface = "bold") +
  labs(y = NULL, subtitle = 'Conditional MM', x="") +
  theme_ggdist() +
  theme(legend.position = "none",
        plot.subtitle = element_text(face="bold", hjust=.5)) + 
  xlim(.25, .65)+
  annotate("label", x = .35, y = 1.5, label = "Mean probability of\nselection among women",
           color="purple", size=2, fontface="bold", hjust=.5)+
  annotate("label", x = .56, y = 1.5, label = "Mean probability of\nselection among men",
           color="#FC4C02", size=2, fontface="bold", hjust=.5)



byGEN.amce<- cj(df_hetero, selected ~ FLAG + BI + TRANS + CONGRUENCE + ETHNICITY + EMOTION +RELATIONSHIP + DISTANCE + EDUCATION_SLIM + JOB_SLIM  + experimental_group,
                  id = ~ ResponseId, estimate="amce",by = ~ resp_gender)

byGEN.amce$lower90 <- byGEN.amce$estimate - (1.645 * byGEN.amce$std.error)
byGEN.amce$higher90 <- byGEN.amce$estimate + (1.645 * byGEN.amce$std.error)
byGEN.amce$lower95 <- byGEN.amce$estimate - (1.96 * byGEN.amce$std.error)
byGEN.amce$higher95 <- byGEN.amce$estimate + (1.96 * byGEN.amce$std.error)
byGEN.amce$lower99 <- byGEN.amce$estimate - (2.58 * byGEN.amce$std.error)
byGEN.amce$higher99 <- byGEN.amce$estimate + (2.58 * byGEN.amce$std.error)
byGEN.amce$lower999 <- byGEN.amce$estimate - (3.29 * byGEN.amce$std.error)
byGEN.amce$higher999 <- byGEN.amce$estimate + (3.29 * byGEN.amce$std.error)

byGEN.amce<-subset(byGEN.amce, feature %in% c("BI")) 
byGEN.amce<-subset(byGEN.amce, level %in% c("Bisexual\nprofile")) 

amceplot2<- ggplot(byGEN.amce, aes(y = level, x = estimate, color = BY)) +
  scale_shape_manual(values = c(10, 18)) +
  scale_fill_manual(values = col2b) +
  scale_color_manual(values = col2b) +
  geom_segment(data = subset(byGEN.amce, BY == "a woman"),
               aes(x = lower90, xend = higher90, y = level, yend = level),
               size = 4, alpha = 1,  position = position_nudge(y = .2)) + 
  geom_segment(data = subset(byGEN.amce, BY == "a woman"),
               aes(x = lower95, xend = higher95, y = level, yend = level),
               size = 4, alpha = .6,  position = position_nudge(y = .2)) + 
  geom_segment(data = subset(byGEN.amce, BY == "a woman"),
               aes(x = lower99, xend = higher99, y = level, yend = level),
               size = 4, alpha = .3,  position = position_nudge(y = .2)) +
  annotate("segment", x = -.051, xend = -.051, y = 1.0, yend = 1.3, color = "grey", linetype = "solid", alpha=.7)+
  annotate("segment", x = -.091, xend = -.051, y = 1.3, yend = 1.3, color = "grey", linetype = "solid", alpha=.7)+
  annotate("segment", x = -.091, xend = -.091, y = 1.2, yend = 1.3, color = "grey", linetype = "solid", alpha=.7)+
  geom_point(data = subset(byGEN.amce, BY == "a woman"),
             aes(x = estimate, y = level), fill = "white", color = "grey67", size = 2.5, pch = 21,  
             position = position_nudge(y = .2))+
  geom_segment(data = subset(byGEN.amce, BY == "a man"),
               aes(x = lower90, xend = higher90, y = level, yend = level),
               size = 4, alpha = 1) + 
  geom_segment(data = subset(byGEN.amce, BY == "a man"),
               aes(x = lower95, xend = higher95, y = level, yend = level),
               size = 4, alpha = .6) + 
  geom_segment(data = subset(byGEN.amce, BY == "a man"),
               aes(x = lower99, xend = higher99, y = level, yend = level),
               size = 4, alpha = .3) +
  geom_point(data = subset(byGEN.amce, BY == "a man"),
             aes(x = estimate, y = level), fill = "white", color = "grey67", size = 2.5, pch = 21) +
  geom_label(data = subset(byGEN.amce, BY == "a woman"),
             aes(label = sprintf("%0.2f", round(estimate, digits = 2)), y = level, x = estimate),
             position = position_nudge(y = .1), size = 2, fontface = "bold") +
  geom_label(data = subset(byGEN.amce, BY == "a man"),
             aes(label = sprintf("%0.2f", round(estimate, digits = 2)), y = level, x = estimate),
             position = position_nudge(y = -.1), size = 2, fontface = "bold") +
  labs(y = NULL, subtitle = 'Conditional AMCE', x="") +
  theme_ggdist() +
  theme(legend.position = "none",
        axis.text.y = element_blank(),
        axis.ticks.y = element_blank(),
        plot.subtitle = element_text(face="bold", hjust=.5)) + 
  xlim(-.12, .02) +
  geom_vline(xintercept = 0, colour = "gray", linetype = "dashed")+
  annotate("label", x = -.07, y = 1.3, label = "0.04 (p<0.01)", size = 2, fontface="bold",color = "black")+
  annotate("label", x = -.08, y = .77, label = "Women",
           color="purple", size=2.5, fontface="bold", hjust=.5)+
  annotate("label", x = -.08, y = .67, label = "Men",
           color="#FC4C02", size=2.5, fontface="bold", hjust=.5)


##Difference in AMCE by gender
coefficients_male <- byGEN.amce$estimate[byGEN.amce$BY == "a man"]
coefficients_female <- byGEN.amce$estimate[byGEN.amce$BY == "a woman"]
SE_male <- byGEN.amce$std.error[byGEN.amce$BY == "a man"]
SE_female <- byGEN.amce$std.error[byGEN.amce$BY == "a woman"]
difference<- coefficients_male-coefficients_female
se_difference <- sqrt(SE_male^2 + SE_female^2)
upper<- difference + (1.96 * se_difference)
lower<- difference - (1.96 * se_difference)
upper99<- difference + (2.58 * se_difference)
lower99<- difference - (2.58 * se_difference)



mmplot2+amceplot2+
  plot_annotation(caption = 'Confidence intervals at 99%, 95% & 90% (N:31,240)',
                  title = 'Study 2: Gender-based results from visual conjoint experiment') & 
  theme(plot.title = element_text(face="bold"))
ggsave(filename = here("figures/Figure8.png"), height=3.5, width=7, dpi = 500)
ggsave(filename = here("figures/tifs/Figure8.tif"), height=3.5, width=7, dpi = 500)




###APPENDIX FIGURES###

bySexuality.mm<- cj(df_full, selected ~ FLAG + BI + TRANS + CONGRUENCE + ETHNICITY + EMOTION +RELATIONSHIP + DISTANCE + EDUCATION_SLIM + JOB_SLIM  + experimental_group,
                    id = ~ ResponseId, estimate="mm",by = ~ queer)

bySexuality.mm$lower90 <- bySexuality.mm$estimate - (1.645 * bySexuality.mm$std.error)
bySexuality.mm$higher90 <- bySexuality.mm$estimate + (1.645 * bySexuality.mm$std.error)
bySexuality.mm$lower95 <- bySexuality.mm$estimate - (1.96 * bySexuality.mm$std.error)
bySexuality.mm$higher95 <- bySexuality.mm$estimate + (1.96 * bySexuality.mm$std.error)
bySexuality.mm$lower99 <- bySexuality.mm$estimate - (2.58 * bySexuality.mm$std.error)
bySexuality.mm$higher99 <- bySexuality.mm$estimate + (2.58 * bySexuality.mm$std.error)
bySexuality.mm$lower999 <- bySexuality.mm$estimate - (3.29 * bySexuality.mm$std.error)
bySexuality.mm$higher999 <- bySexuality.mm$estimate + (3.29 * bySexuality.mm$std.error)

bySexuality.mm<-subset(bySexuality.mm, feature %in% c("BI")) 
bySexuality.mm <- bySexuality.mm %>%
  mutate(level = ifelse(level == "Heterosexual\nprofile", "Non-bisexual\nprofile", "Bisexual\nprofile"))

col2 <- c("steelblue","#D2042D")

bySex_mmplot<- ggplot(bySexuality.mm, aes(y = level, x = estimate, color = BY)) +
  scale_shape_manual(values = c(10, 18)) +
  scale_fill_manual(values = col2) +
  scale_color_manual(values = col2) +
  geom_vline(xintercept = 0.45, colour="gray", linetype="dashed")+ ##.45 median acceptance
  geom_segment(data = subset(bySexuality.mm, BY == "Heterosexual"),
               aes(x = lower90, xend = higher90, y = level, yend = level),
               size = 4, alpha = 1,  position = position_nudge(y = .2)) + 
  geom_segment(data = subset(bySexuality.mm, BY == "Heterosexual"),
               aes(x = lower95, xend = higher95, y = level, yend = level),
               size = 3, alpha = .8,  position = position_nudge(y = .2)) + 
  geom_segment(data = subset(bySexuality.mm, BY == "Heterosexual"),
               aes(x = lower99, xend = higher99, y = level, yend = level),
               size = 2, alpha = .6,  position = position_nudge(y = .2)) +
  geom_point(data = subset(bySexuality.mm, BY == "Heterosexual"),
             aes(x = estimate, y = level), fill = "white", color = "grey67", size = 2.5, pch = 21,  
             position = position_nudge(y = .2))+
  geom_segment(data = subset(bySexuality.mm, BY == "LGBTQ+"),
               aes(x = lower90, xend = higher90, y = level, yend = level),
               size = 4, alpha = 1) + 
  geom_segment(data = subset(bySexuality.mm, BY == "LGBTQ+"),
               aes(x = lower95, xend = higher95, y = level, yend = level),
               size = 3, alpha = .8) + 
  geom_segment(data = subset(bySexuality.mm, BY == "LGBTQ+"),
               aes(x = lower99, xend = higher99, y = level, yend = level),
               size = 2, alpha = .6) +
  geom_point(data = subset(bySexuality.mm, BY == "LGBTQ+"),
             aes(x = estimate, y = level), fill = "white", color = "grey67", size = 2.5, pch = 21) +
  geom_label(data = subset(bySexuality.mm, BY == "Heterosexual"),
             aes(label = sprintf("%0.2f", round(estimate, digits = 2)), y = level, x = estimate),
             position = position_nudge(y = .1), size = 3, fontface = "bold"
  ) +
  geom_label(data = subset(bySexuality.mm, BY == "LGBTQ+"),
             aes(label = sprintf("%0.2f", round(estimate, digits = 2)), y = level, x = estimate),
             position = position_nudge(y = -.1), size = 3, fontface = "bold") +
  labs(y = NULL, subtitle = 'Conditional MM by sexual & gender identity', x="") +
  theme_ggdist() +
  theme(legend.position = "none",
        plot.subtitle = element_text(face="bold", hjust=.5))+
  annotate("label", x = .4, y = 2, label = "Heterosexual respondent",
           color="steelblue", size=3.5, fontface="bold", hjust=.5)+
  annotate("label", x = .4, y = 2.1, label = "LGBTQ+ respondent",
           color="#D2042D", size=3.5, fontface="bold", hjust=.5)
ggsave( filename = here("figures/appendix/FigureB1.png"), dpi = 300)



df_full <- df_full %>%
  mutate(combined = case_when(
    gender == "a man" & queer == "Heterosexual" ~ "Heterosexual man",
    gender == "a man" & queer == "LGBTQ+" ~ "LGBTQ+ man",
    gender == "a woman" & queer == "Heterosexual" ~ "Heterosexual woman",
    gender == "a woman" & queer == "LGBTQ+" ~ "LGBTQ+ woman"
  ))
df_full <- df_full %>%
  mutate(combined =as.factor(combined))


by4arm.mm<- cj(df_full, selected ~ FLAG + BI + TRANS + CONGRUENCE + ETHNICITY + EMOTION +RELATIONSHIP + DISTANCE + EDUCATION_SLIM + JOB_SLIM + experimental_group,
               id = ~ ResponseId, estimate="mm",by = ~ combined)
by4arm.mm<-subset(by4arm.mm, feature %in% c("BI")) 

bySexGender_mmplot <- plot(by4arm.mm, feature_headers = FALSE) +
  ggplot2::aes(x = estimate, y = level, colour = BY, shape=BY) + 
  ggplot2::facet_wrap(~feature, ncol = 1L, labeller=labeller(feature = labels),
                      scales = "free_y", strip.position = "right")+
  labs(x="",  subtitle="Conditional MM")+
  geom_vline(xintercept = 0.35, colour = "grey", linetype = "dashed")+
  geom_vline(xintercept = 0.56, colour = "grey", linetype = "dashed")+
  scale_color_manual(values = c("blue1", "blue4", "red1", "red4"))+
  theme_ggdist()+
  guides(col = guide_legend(nrow = 2))+
  theme(strip.text = element_blank(),
        legend.title = element_blank(),
        legend.position="bottom",
        legend.text = element_text(face="bold"),
        axis.title.y = element_blank(),
        axis.text.y = element_blank(),
        axis.ticks.y = element_blank(),
        plot.subtitle = element_text(face="bold", hjust=.5))


by4arm.amce<- cj(df_full, selected ~ FLAG + BI + TRANS + CONGRUENCE + ETHNICITY + EMOTION +RELATIONSHIP + DISTANCE + EDUCATION_SLIM + JOB_SLIM  + experimental_group,
                 id = ~ ResponseId, estimate="amce",by = ~ combined)
by4arm.amce<-subset(by4arm.amce, feature %in% c("BI")) 

bySexGender_amceplot <- plot(by4arm.amce, feature_headers = FALSE) +
  ggplot2::aes(x = estimate, y = level, colour = BY, shape=BY) + 
  ggplot2::facet_wrap(~feature, ncol = 1L, labeller=labeller(feature = labels),
                      scales = "free_y", strip.position = "right")+
  labs(x="",  subtitle="Conditional AMCE")+
  geom_vline(xintercept = 0, colour = "grey", linetype = "dashed")+
  scale_color_manual(values = c("blue1", "blue4", "red1", "red4"))+
  theme_ggdist()+
  guides(col = guide_legend(nrow = 2))+
  theme(strip.text = element_blank(),
        legend.title = element_blank(),
        legend.position="bottom",
        legend.text = element_text(face="bold"),
        axis.title.y = element_blank(),
        axis.text.y = element_blank(),
        axis.ticks.y = element_blank(),
        plot.subtitle = element_text(face="bold", hjust=.5))
bySexGender_mmplot+bySexGender_amceplot
ggsave("figures/appendix/FigureB2.png", dpi = 600)


fourway.mm<- cj(df_hetero, selected ~ FLAG + BI + TRANS + CONGRUENCE + ETHNICITY + EMOTION +RELATIONSHIP + DISTANCE + EDUCATION_SLIM + JOB_SLIM  + experimental_group,
                id = ~ ResponseId, estimate="mm",by = ~ derecha_mujer)

fourway.mm$lower90 <- fourway.mm$estimate - (1.645 * fourway.mm$std.error)
fourway.mm$higher90 <- fourway.mm$estimate + (1.645 * fourway.mm$std.error)
fourway.mm$lower95 <- fourway.mm$estimate - (1.96 * fourway.mm$std.error)
fourway.mm$higher95 <- fourway.mm$estimate + (1.96 * fourway.mm$std.error)
fourway.mm$lower99 <- fourway.mm$estimate - (2.58 * fourway.mm$std.error)
fourway.mm$higher99 <- fourway.mm$estimate + (2.58 * fourway.mm$std.error)
fourway.mm$lower999 <- fourway.mm$estimate - (3.29 * fourway.mm$std.error)
fourway.mm$higher999 <- fourway.mm$estimate + (3.29 * fourway.mm$std.error)

fourway.mm<-subset(fourway.mm, feature %in% c("BI")) 

MM_4way<- ggplot(fourway.mm, aes(y = level, x = estimate, color = BY)) +
  scale_shape_manual(values = c(10, 18)) +
  scale_fill_manual(values = col4) +
  scale_color_manual(values = col4)+ 
  geom_vline(xintercept = 0.35, colour = "purple", linetype = "dashed")+
  geom_vline(xintercept = 0.56, colour = "#FC4C02", linetype = "dashed")+
  geom_segment(data = subset(fourway.mm, BY == "Left-wing woman"),
               aes(x = lower90, xend = higher90, y = level, yend = level),
               size = 4, alpha = 1) + 
  geom_segment(data = subset(fourway.mm, BY == "Left-wing woman"),
               aes(x = lower95, xend = higher95, y = level, yend = level),
               size = 3, alpha = .8) + 
  geom_segment(data = subset(fourway.mm, BY == "Left-wing woman"),
               aes(x = lower99, xend = higher99, y = level, yend = level),
               size = 2, alpha = .6) +
  geom_point(data = subset(fourway.mm, BY == "Left-wing woman"),
             aes(x = estimate, y = level), fill = "white", color = "grey67", size = 2.5, pch = 21)+
  geom_segment(data = subset(fourway.mm, BY == "Right-wing woman"),
               aes(x = lower90, xend = higher90, y = level, yend = level),
               size = 4, alpha = 1,  position = position_nudge(y = .4)) + 
  geom_segment(data = subset(fourway.mm, BY == "Right-wing woman"),
               aes(x = lower95, xend = higher95, y = level, yend = level),
               size = 3, alpha = .8,  position = position_nudge(y = .4)) + 
  geom_segment(data = subset(fourway.mm, BY == "Right-wing woman"),
               aes(x = lower99, xend = higher99, y = level, yend = level),
               size = 2, alpha = .6,  position = position_nudge(y = .4)) +
  geom_point(data = subset(fourway.mm, BY == "Right-wing woman"),
             aes(x = estimate, y = level), fill = "white", color = "grey67", size = 2.5, pch = 21,  
             position = position_nudge(y = .4))+
  geom_segment(data = subset(fourway.mm, BY == "Left-wing man"),
               aes(x = lower90, xend = higher90, y = level, yend = level),
               size = 4, alpha = 1) + 
  geom_segment(data = subset(fourway.mm, BY == "Left-wing man"),
               aes(x = lower95, xend = higher95, y = level, yend = level),
               size = 3, alpha = .8) + 
  geom_segment(data = subset(fourway.mm, BY == "Left-wing man"),
               aes(x = lower99, xend = higher99, y = level, yend = level),
               size = 2, alpha = .6) +
  geom_point(data = subset(fourway.mm, BY == "Left-wing man"),
             aes(x = estimate, y = level), fill = "white", color = "grey67", size = 2.5, pch = 21) +
  geom_segment(data = subset(fourway.mm, BY == "Right-wing man"),
               aes(x = lower90, xend = higher90, y = level, yend = level),
               size = 4, alpha = 1,  position = position_nudge(y = .4)) + 
  geom_segment(data = subset(fourway.mm, BY == "Right-wing man"),
               aes(x = lower95, xend = higher95, y = level, yend = level),
               size = 3, alpha = .8,  position = position_nudge(y = .4)) + 
  geom_segment(data = subset(fourway.mm, BY == "Right-wing man"),
               aes(x = lower99, xend = higher99, y = level, yend = level),
               size = 2, alpha = .6,  position = position_nudge(y = .4)) +
  geom_point(data = subset(fourway.mm, BY == "Right-wing man"),
             aes(x = estimate, y = level), fill = "white", color = "grey67", size = 2.5, pch = 21,  
             position = position_nudge(y = .4))+
  labs(y = NULL, subtitle = 'Conditional MM', x="") +
  theme_ggdist() +
  geom_label(data = subset(fourway.mm, BY == "Left-wing woman"),
             aes(label = sprintf("%0.2f", round(estimate, digits = 2)), y = level, x = estimate),
             position = position_nudge(y = -.1), size = 3, fontface = "bold", show.legend = FALSE) +
  geom_label(data = subset(fourway.mm, BY == "Right-wing woman"),
             aes(label = sprintf("%0.2f", round(estimate, digits = 2)), y = level, x = estimate),
             position = position_nudge(y = .3), size = 3, fontface = "bold", show.legend = FALSE) +
  geom_label(data = subset(fourway.mm, BY == "Left-wing man"),
             aes(label = sprintf("%0.2f", round(estimate, digits = 2)), y = level, x = estimate),
             position = position_nudge(y = -.1), size = 3, fontface = "bold", show.legend = FALSE) +
  geom_label(data = subset(fourway.mm, BY == "Right-wing man"),
             aes(label = sprintf("%0.2f", round(estimate, digits = 2)), y = level, x = estimate),
             position = position_nudge(y = .3), size = 3, fontface = "bold", show.legend = FALSE)+
  theme(legend.position = "bottom",
        legend.title=element_blank(),
        plot.subtitle = element_text(face="bold", hjust=.5))+
  guides(color = guide_legend(nrow = 2))


fourway.amce<- cj(df_hetero, selected ~ FLAG + BI + TRANS + CONGRUENCE + ETHNICITY + EMOTION +RELATIONSHIP + DISTANCE + EDUCATION_SLIM + JOB_SLIM  + experimental_group,
                  id = ~ ResponseId, estimate="amce",by = ~ derecha_mujer)

fourway.amce$lower90 <- fourway.amce$estimate - (1.645 * fourway.amce$std.error)
fourway.amce$higher90 <- fourway.amce$estimate + (1.645 * fourway.amce$std.error)
fourway.amce$lower95 <- fourway.amce$estimate - (1.96 * fourway.amce$std.error)
fourway.amce$higher95 <- fourway.amce$estimate + (1.96 * fourway.amce$std.error)
fourway.amce$lower99 <- fourway.amce$estimate - (2.58 * fourway.amce$std.error)
fourway.amce$higher99 <- fourway.amce$estimate + (2.58 * fourway.amce$std.error)
fourway.amce$lower999 <- fourway.amce$estimate - (3.29 * fourway.amce$std.error)
fourway.amce$higher999 <- fourway.amce$estimate + (3.29 * fourway.amce$std.error)

fourway.amce<-subset(fourway.amce, feature %in% c("BI")) 
fourway.amce<-subset(fourway.amce, level %in% c("Bisexual\nprofile")) 


AMCE_4way<- ggplot(fourway.amce, aes(y = level, x = estimate, color = BY)) +
  scale_fill_manual(values = col4) +
  scale_color_manual(values = col4) +
  geom_segment(data = subset(fourway.amce, BY == "Left-wing woman"),
               aes(x = lower90, xend = higher90, y = level, yend = level),
               size = 4, alpha = 1) + 
  geom_segment(data = subset(fourway.amce, BY == "Left-wing woman"),
               aes(x = lower95, xend = higher95, y = level, yend = level),
               size = 3, alpha = .8) + 
  geom_segment(data = subset(fourway.amce, BY == "Left-wing woman"),
               aes(x = lower99, xend = higher99, y = level, yend = level),
               size = 2, alpha = .6) +
  geom_point(data = subset(fourway.amce, BY == "Left-wing woman"),
             aes(x = estimate, y = level), fill = "white", color = "grey67", size = 2.5, pch = 21)+
  geom_segment(data = subset(fourway.amce, BY == "Right-wing woman"),
               aes(x = lower90, xend = higher90, y = level, yend = level),
               size = 4, alpha = 1,  position = position_nudge(y = .2)) + 
  geom_segment(data = subset(fourway.amce, BY == "Right-wing woman"),
               aes(x = lower95, xend = higher95, y = level, yend = level),
               size = 3, alpha = .8,  position = position_nudge(y = .2)) + 
  geom_segment(data = subset(fourway.amce, BY == "Right-wing woman"),
               aes(x = lower99, xend = higher99, y = level, yend = level),
               size = 2, alpha = .6,  position = position_nudge(y = .2)) +
  geom_point(data = subset(fourway.amce, BY == "Right-wing woman"),
             aes(x = estimate, y = level), fill = "white", color = "grey67", size = 2.5, pch = 21,  
             position = position_nudge(y = .2))+
  geom_segment(data = subset(fourway.amce, BY == "Left-wing man"),
               aes(x = lower90, xend = higher90, y = level, yend = level),
               size = 4, alpha = 1, position = position_nudge(y = -.2)) + 
  geom_segment(data = subset(fourway.amce, BY == "Left-wing man"),
               aes(x = lower95, xend = higher95, y = level, yend = level),
               size = 3, alpha = .8, position = position_nudge(y = -.2)) + 
  geom_segment(data = subset(fourway.amce, BY == "Left-wing man"),
               aes(x = lower99, xend = higher99, y = level, yend = level),
               size = 2, alpha = .6, position = position_nudge(y = -.2)) +
  geom_point(data = subset(fourway.amce, BY == "Left-wing man"),
             aes(x = estimate, y = level), fill = "white", color = "grey67", size = 2.5, pch = 21,
             position = position_nudge(y = -.2)) +
  geom_segment(data = subset(fourway.amce, BY == "Right-wing man"),
               aes(x = lower90, xend = higher90, y = level, yend = level),
               size = 4, alpha = 1,  position = position_nudge(y = -.4)) + 
  geom_segment(data = subset(fourway.amce, BY == "Right-wing man"),
               aes(x = lower95, xend = higher95, y = level, yend = level),
               size = 3, alpha = .8,  position = position_nudge(y = -.4)) + 
  geom_segment(data = subset(fourway.amce, BY == "Right-wing man"),
               aes(x = lower99, xend = higher99, y = level, yend = level),
               size = 2, alpha = .6,  position = position_nudge(y = -.4)) +
  geom_point(data = subset(fourway.amce, BY == "Right-wing man"),
             aes(x = estimate, y = level), fill = "white", color = "grey67", size = 2.5, pch = 21,  
             position = position_nudge(y = -.4))+
  theme_ggdist() +
  geom_label(data = subset(fourway.amce, BY == "Left-wing woman"),
             aes(label = sprintf("%0.2f", round(estimate, digits = 2)), y = level, x = estimate),
             position = position_nudge(y = -.1), size = 3, fontface = "bold", show.legend = FALSE) +
  geom_label(data = subset(fourway.amce, BY == "Right-wing woman"),
             aes(label = sprintf("%0.2f", round(estimate, digits = 2)), y = level, x = estimate),
             position = position_nudge(y = .1), size = 3, fontface = "bold", show.legend = FALSE) +
  geom_label(data = subset(fourway.amce, BY == "Left-wing man"),
             aes(label = sprintf("%0.2f", round(estimate, digits = 2)), y = level, x = estimate),
             position = position_nudge(y = -.3), size = 3, fontface = "bold", show.legend = FALSE) +
  geom_label(data = subset(fourway.amce, BY == "Right-wing man"),
             aes(label = sprintf("%0.2f", round(estimate, digits = 2)), y = level, x = estimate),
             position = position_nudge(y = -.5), size = 3, fontface = "bold", show.legend = FALSE)+
  labs(x="", y="",
       subtitle="Conditional AMCE")+
  theme(legend.position = "bottom",
        legend.title=element_blank(),
        axis.text.y = element_blank(),
        axis.ticks.y = element_blank(),
        plot.subtitle = element_text(face="bold", hjust=.5)) + 
  guides(color = guide_legend(nrow = 2))+
  geom_vline(xintercept = 0, colour = "gray", linetype = "dashed")


MM_4way+AMCE_4way+
  plot_annotation(caption = 'Confidence intervals at 99%, 95% & 90% (N:31,240)',
                  title = 'Study 2: Four-way subgroup analysis') & 
  theme(plot.title = element_text(face="bold"))
ggsave(filename = here("figures/appendix/FigureB3.png"),height=6, width=9, dpi = 300)
